home *** CD-ROM | disk | FTP | other *** search
- #ifndef lint
- static char *RCSid = "$Id: hidden3d.c,v 1.6 1995/12/02 22:04:34 drd Exp $";
- #endif
-
-
- /* GNUPLOT - hidden3d.c */
- /*
- * Copyright (C) 1986 - 1993 Thomas Williams, Colin Kelley
- *
- * Permission to use, copy, and distribute this software and its
- * documentation for any purpose with or without fee is hereby granted,
- * provided that the above copyright notice appear in all copies and
- * that both that copyright notice and this permission notice appear
- * in supporting documentation.
- *
- * Permission to modify the software is granted, but not the right to
- * distribute the modified code. Modifications are to be distributed
- * as patches to released version.
- *
- * This software is provided "as is" without express or implied warranty.
- *
- *
- * AUTHORS
- *
- * Original Software:
- * Gershon Elber and many others.
- *
- * 19 September 1992 Lawrence Crowl (crowl@cs.orst.edu)
- * Added user-specified bases for log scaling.
- *
- * 3.6 - split graph3d.c into graph3d.c (graph),
- * util3d.c (intersections, etc)
- * hidden3d.c (hidden-line removal code)
- *
- * There is a mailing list for gnuplot users. Note, however, that the
- * newsgroup
- * comp.graphics.gnuplot
- * is identical to the mailing list (they
- * both carry the same set of messages). We prefer that you read the
- * messages through that newsgroup, to subscribing to the mailing list.
- * (If you can read that newsgroup, and are already on the mailing list,
- * please send a message info-gnuplot-request@dartmouth.edu, asking to be
- * removed from the mailing list.)
- *
- * The address for mailing to list members is
- * info-gnuplot@dartmouth.edu
- * and for mailing administrative requests is
- * info-gnuplot-request@dartmouth.edu
- * The mailing list for bug reports is
- * bug-gnuplot@dartmouth.edu
- * The list of those interested in beta-test versions is
- * info-gnuplot-beta@dartmouth.edu
- */
-
- #include <math.h>
- #if !defined(sequent) && !defined(apollo) && !defined(alliant)
- #include <limits.h>
- #endif
- #include "plot.h"
- #include "setshow.h"
-
- extern int suppressMove;
-
- extern int xright,xleft,ybot,ytop;
-
- extern double min_array[], max_array[];
- extern int auto_array[], log_array[];
- extern double base_array[], log_base_array[];
-
- /* for convenience while converting to use these arrays */
- #define x_min3d min_array[FIRST_X_AXIS]
- #define x_max3d max_array[FIRST_X_AXIS]
- #define y_min3d min_array[FIRST_Y_AXIS]
- #define y_max3d max_array[FIRST_Y_AXIS]
- #define z_min3d min_array[FIRST_Z_AXIS]
- #define z_max3d max_array[FIRST_Z_AXIS]
- #define min3d_z min_array[FIRST_Z_AXIS]
- #define max3d_z max_array[FIRST_Z_AXIS]
-
- extern int base_z;
- extern int hidden_no_update, hidden_active;
- extern int hidden_line_type_above, hidden_line_type_below;
-
- static int zsort __P((int * r1, int * r2));
-
-
- /* We divvy up the figure into the component boxes that make it up, and then
- sort them by the z-value (which is really just an average value). */
- struct pnts{
- unsigned int x,y,z;
- int flag;
- long int style_used; /* acw test */
- int nplot;
- };
- static int * boxlist;
- static struct pnts * nodes;
-
- /* These variables are used to keep track of the range of x values used in the
- line drawing routine. */
-
- /* THESE SHOULD BE STATIC but util3d.c needs them - try to move code
- * from util3d.c into here
- */
-
- long int xmin_hl,xmax_hl;
- /* These arrays are used to keep track of the minimum and maximum y values used
- for each X value. These are only used for drawing the individual boxes that
- make up the 3d figure. After each box is drawn, the information is copied
- to the bitmap. */
- short int *ymin_hl=0, *ymax_hl=0;
- /*
- * These numbers are chosen as dividers into the bitmap.
- */
- int xfact, yfact;
- #define XREDUCE(X) ((X)/xfact)
- #define YREDUCE(Y) ((Y)/yfact)
- /* Bitmap of the screen. The array for each x value is malloc-ed as needed */
- unsigned short int **pnt=0;
- #define IFSET(X,Y) (pnt[X] == 0 ? 0 : (((pnt[X])[(Y)>>4] >> ((Y) & 0xf)) & 0x01))
-
-
- /* Initialize the necessary steps for hidden line removal. */
- void init_hidden_line_removal()
- {
- int i;
- /* We want to keep the bitmap size less than 2048x2048, so we choose
- * integer dividers for the x and y coordinates to keep the x and y
- * ranges less than 2048. In practice, the x and y sizes for the bitmap
- * will be somewhere between 1024 and 2048, except in cases where the
- * coordinates ranges for the device are already less than 1024.
- * We do this mainly to control the size of the bitmap, but it also
- * speeds up the computation. We maintain separate dividers for
- * x and y.
- */
- xfact = (xright-xleft)/1024;
- yfact = (ytop-ybot)/1024;
- if(xfact == 0) xfact=1;
- if(yfact == 0) yfact=1;
- if(pnt == 0){
- i = sizeof(short int*)*(XREDUCE(xright) - XREDUCE(xleft) + 1);
- pnt = (unsigned short int **) alloc((unsigned long)i, "hidden");
- memset(pnt, 0, i);
- };
- ymin_hl = (short int *) alloc((unsigned long)sizeof(short int)*
- (XREDUCE(xright) - XREDUCE(xleft) + 1), "hidden");
- ymax_hl = (short int *) alloc((unsigned long)sizeof(short int)*
- (XREDUCE(xright) - XREDUCE(xleft) + 1), "hidden");
- }
-
- /* Reset the hidden line data to a fresh start. */
- void reset_hidden_line_removal()
- {
- int i;
- if(pnt){
- for(i=0;i<=XREDUCE(xright)-XREDUCE(xleft);i++) {
- if(pnt[i])
- { free(pnt[i]); pnt[i] = 0;};
- };
- };
- }
-
- /* Terminates the hidden line removal process. Free any memory allocated by */
- /* init_hidden_line_removal above. */
- void term_hidden_line_removal()
- {
- if(pnt){
- int j;
- for(j=0;j<=XREDUCE(xright)-XREDUCE(xleft);j++) {
- if(pnt[j])
- { free(pnt[j]); pnt[j] = 0;};
- };
- free(pnt);
- pnt = 0;
- };
- if(ymin_hl) free(ymin_hl), ymin_hl = 0;
- if(ymax_hl) free(ymax_hl), ymax_hl = 0;
- }
-
-
- static int zsort( r1, r2)
- int * r1;
- int * r2;
- {
- int z1, z2;
- z1 = nodes[*r1].z;
- z2 = nodes[*r2].z;
- if (z1 < z2) return 1;
- if (z1 == z2) return 0;
- return -1;
- }
- #define TESTBOX(X,Y) \
- if(X<xmin_box) xmin_box = X; \
- if(X>xmax_box) xmax_box = X; \
- if(Y<ymin_box) ymin_box = Y; \
- if(Y>ymax_box) ymax_box = Y;
- /* Usefull macro to help us figure out which side of the surface we are on */
- #define XPRD(I,J,K) \
- ((nodes[I].x-nodes[J].x)*(nodes[J].y-nodes[K].y) - \
- (nodes[I].y-nodes[J].y)*(nodes[J].x-nodes[K].x))
- #define MAYBE_LINEPOINT(J) \
- if((nodes[J].flag & 0x20) != 0) { \
- x = nodes[J].x; \
- y = nodes[J].y; \
- nodes[J].flag -= 0x20; \
- if (!clip_point(x,y) && \
- !IFSET(XREDUCE(x)-XREDUCE(xleft),YREDUCE(y)-YREDUCE(ybot))) \
- (*t->point)(x,y, plot_info[nplot].point_type); \
- };
-
- struct surface_plots{
- int above_color;
- int below_color;
- int row_offset;
- int point_type;
- };
- /* All of the plots coming into this routine are assumed to have grid
- topology. */
-
- void plot3d_hidden(plots, pcount)
- struct surface_points *plots;
- int pcount;
- {
- struct surface_points *this_plot;
- long int i, j;
- int nplot;
- long int x,y,z ,nseg, ncrv, ncrv1; /* point in terminal coordinates */
- unsigned short int * cpnt;
- short int mask1, mask2;
- long int indx1, indx2, k, m;
- short int xmin_box, xmax_box, ymin_box, ymax_box;
- struct surface_plots * plot_info;
- int row_offset, nnode;
- short int y_malloc; /* Amount of space we need for one vertical row of
- bitmap, and byte offset of first used element */
- struct termentry *t = term;
- struct iso_curve *icrvs;
- int current_style = 0x7fff; /* Current line style */
- int surface;
- nnode = 0;
- nseg = 0;
- nplot = 0;
- this_plot = plots;
-
- for (surface = 0;
- surface < pcount;
- this_plot = this_plot->next_sp, surface++) {
- nplot++;
- icrvs = plots->iso_crvs;
- icrvs = plots->iso_crvs;
- if(this_plot->plot_type == FUNC3D) {
- ncrv=0;
- for(icrvs = this_plot->iso_crvs;icrvs;icrvs=icrvs->next)
- ++ncrv;
- }
- else if(this_plot->plot_type == DATA3D)
- ncrv = this_plot->num_iso_read;
- else {
- graph_error("Plot type is neither function nor data");
- return; /* stops a gcc -Wunitialised warning */
- }
- nnode += ncrv * (this_plot->iso_crvs->p_count);
- /* for(icrvs = this_plot->iso_crvs,ncrv=0;icrvs;icrvs=icrvs->next,ncrv++) { };
- nnode += ncrv * (this_plot->iso_crvs->p_count); */
- switch(this_plot->plot_style) {
- case YERRORBARS:
- case XERRORBARS:
- case XYERRORBARS:
- case DOTS:
- case POINTSTYLE:
- case LINESPOINTS:
- nseg += (ncrv) * (this_plot->iso_crvs->p_count);
- break;
- case LINES:
- nseg += (ncrv-1) * (this_plot->iso_crvs->p_count-1);
- break;
- case IMPULSES:
- /* There will be two nodes for each segment */
- nnode += ncrv * (this_plot->iso_crvs->p_count);
- nseg += (ncrv) * (this_plot->iso_crvs->p_count);
- break;
- }
- };
- boxlist = (int *) alloc((unsigned long)sizeof(int)*nseg, "hidden");
- nodes = (struct pnts *) alloc((unsigned long)sizeof(struct pnts)*nnode, "hidden");
- plot_info = (struct surface_plots *) alloc((unsigned long)sizeof(struct surface_plots)*nplot,"hidden");
- nnode = 0;
- nseg = 0;
- nplot = 0;
- this_plot = plots;
- hidden_no_update = FALSE;
-
- if ( hidden3d && draw_surface)
- for (surface = 0;
- surface < pcount;
- this_plot = this_plot->next_sp, surface++) {
- (*t->linetype)(this_plot->line_type);
- hidden_line_type_above = this_plot->line_type;
- hidden_line_type_below = this_plot->line_type + 1;
- if(this_plot->plot_type == FUNC3D) {
- for(icrvs = this_plot->iso_crvs,ncrv=0;icrvs;icrvs=icrvs->next,ncrv++) { };
- /* if(this_plot->has_grid_topology) ncrv >>= 1; */
- };
- if(this_plot->plot_type == DATA3D)
- ncrv = this_plot->num_iso_read;
- icrvs = this_plot->iso_crvs;
- ncrv1 = ncrv;
- ncrv = 0;
- while ( icrvs) {
- struct coordinate GPHUGE *points = icrvs->points;
- for (i = 0; i < icrvs->p_count; i++) {
- map3d_xy(points[i].x, points[i].y, points[i].z,&nodes[nnode].x,&nodes[nnode].y);
- nodes[nnode].z = map3d_z(points[i].x, points[i].y, points[i].z);
- nodes[nnode].flag = (i==0 ? 1 : 0) + (ncrv == 0 ? 2 : 0) +
- (i == icrvs->p_count-1 ? 4 : 0) + (ncrv == ncrv1-1 ? 8 : 0);
- nodes[nnode].nplot = nplot;
- nodes[nnode].style_used = -1000; /* indicates no style */
- switch(this_plot->plot_style) {
- case LINESPOINTS:
- if(i < icrvs->p_count-1 && ncrv < ncrv1-1)
- nodes[nnode].flag |= 0x30;
- else
- nodes[nnode].flag |= 0x20;
- boxlist[nseg++] = nnode++;
- break;
- case LINES:
- if(i < icrvs->p_count-1 && ncrv < ncrv1-1)
- {
- nodes[nnode].flag |= 0x10;
- boxlist[nseg++] = nnode++;
- }
- else
- nnode++;
- break;
- case YERRORBARS:
- case XERRORBARS:
- case XYERRORBARS:
- case POINTSTYLE:
- case DOTS:
- nodes[nnode].flag |= 0x40;
- boxlist[nseg++] = nnode++;
- break;
- case IMPULSES:
- nodes[nnode].flag |= 0x80;
- boxlist[nseg++] = nnode++;
- map3d_xy(points[i].x, points[i].y, (double)base_z,
- &nodes[nnode].x,&nodes[nnode].y);
- nodes[nnode].z = map3d_z(points[i].x, points[i].y, (double)base_z);
- nnode++;
- break;
- }
- }
- icrvs = icrvs->next;
- ncrv++;
- if(ncrv == ncrv1) break;
- }
- /* Next we go through all of the boxes, and substitute the average z value
- for the box for the z value of the corner node */
- plot_info[nplot].above_color = this_plot->line_type;
- plot_info[nplot].below_color = this_plot->line_type+1;
- plot_info[nplot].point_type =
- ((this_plot->plot_style == DOTS) ? -1 : this_plot->point_type);
- plot_info[nplot++].row_offset = this_plot->iso_crvs->p_count;
- }
- for(i=0; i<nseg; i++){
- j = boxlist[i];
- if ((nodes[j].flag & 0x80) != 0) {
- nodes[j].z = (nodes[j].z < nodes[j+1].z ? nodes[j].z : nodes[j+1].z);
- continue;
- };
- if ((nodes[j].flag & 0x10) == 0) continue;
- row_offset = plot_info[nodes[j].nplot].row_offset;
- z = nodes[j].z;
- if (z < nodes[j+1].z) z = nodes[j+1].z;
- if (z < nodes[j+row_offset].z) z = nodes[j+row_offset].z;
- if (z < nodes[j+row_offset+1].z) z = nodes[j+row_offset+1].z;
- };
- qsort (boxlist, nseg, sizeof(int), (sortfunc)zsort);
- y_malloc = (2+ (YREDUCE(ytop)>>4) - (YREDUCE(ybot)>>4))*sizeof(short int);
- for(i=0;i<=(XREDUCE(xright)-XREDUCE(xleft));i++) {
- ymin_hl[i] = 0x7fff;
- ymax_hl[i] = 0;
- };
- for(i=0;i<nseg;i++) {
- j = boxlist[i];
- nplot = nodes[j].nplot;
- row_offset = plot_info[nplot].row_offset;
- if((nodes[j].flag & 0x40) != 0) {
- x = nodes[j].x;
- y = nodes[j].y;
- if (!clip_point(x,y) &&
- !IFSET(XREDUCE(x)-XREDUCE(xleft),YREDUCE(y)-YREDUCE(ybot)))
- (*t->point)(x,y, plot_info[nplot].point_type);
- };
- if((nodes[j].flag & 0x80) != 0) { /* impulses */
- clip_move(nodes[j].x,nodes[j].y);
- clip_vector(nodes[j+1].x,nodes[j+1].y);
- };
- if((nodes[j].flag & 0x10) != 0) {
- /* It is possible, and often profitable, to take a quick look and see
- if the current box is entirely obscured. If this is the case we will
- not even bother testing this box any further. */
- xmin_box = 0x7fff;
- xmax_box = 0;
- ymin_box = 0x7fff;
- ymax_box = 0;
- TESTBOX(nodes[j].x-xleft,nodes[j].y-ybot);
- TESTBOX(nodes[j+1].x-xleft,nodes[j+1].y-ybot);
- TESTBOX(nodes[j+row_offset].x-xleft,nodes[j+row_offset].y-ybot);
- TESTBOX(nodes[j+row_offset+1].x-xleft,nodes[j+row_offset+1].y-ybot);
- z=0;
- if(xmin_box < 0) xmin_box = 0;
- if(ymin_box < 0) ymin_box = 0;
- if(xmax_box > xright-xleft) xmax_box = xright-xleft;
- if(ymax_box > ytop-ybot) ymax_box = ytop-ybot;
- /* Now check bitmap. These coordinates have not been reduced */
- if(xmin_box <= xmax_box && ymin_box <= ymax_box){
- ymin_box = YREDUCE(ymin_box);
- ymax_box = YREDUCE(ymax_box);
- xmin_box = XREDUCE(xmin_box);
- xmax_box = XREDUCE(xmax_box);
- indx1 = ymin_box >> 4;
- indx2 = ymax_box >> 4;
- mask1 = 0xffff << (ymin_box & 0x0f);
- mask2 = 0xffff >> (0x0f-(ymax_box & 0x0f));
- for(m=xmin_box;m<=xmax_box;m++) {
- if(pnt[m] == 0) {z++; break;};
- cpnt = pnt[m] + indx1;
- if(indx1 == indx2){
- if((*cpnt & mask1 & mask2) != (mask1 & mask2)) {z++; break;}
- } else {
- if((*cpnt++ & mask1) != mask1) {z++; break;}
- k = indx1+1;
- while (k != indx2) {
- if(*cpnt++ != 0xffff) {z++; break;}
- k++;
- };
- if((*cpnt++ & mask2) != mask2) {z++; break;}
- };
- };
- };
- /* z is 0 if all of the pixels used by the current box are already covered.
- No point in proceeding, so we just skip all further processing of this
- box. */
- if(!z) continue;
- /* Now we need to figure out whether we are looking at the top or the
- bottom of the square. A simple cross product will tell us this.
- If the square is really distorted then this will not be accurate,
- but in such cases we would actually be seeing both sides at the same
- time. We choose the vertex with the largest z component to
- take the cross product at. */
- {
- int z1, z2 ,z3, z4;
- z1 = XPRD(j+row_offset,j,j+1);
- z2 = XPRD(j,j+1,j+1+row_offset);
- z3 = XPRD(j+1,j+row_offset+1,j+row_offset);
- z4 = XPRD(j+row_offset+1,j+row_offset,j);
- z=0;
- z += (z1 > 0 ? 1 : -1);
- z += (z2 > 0 ? 1 : -1);
- z += (z3 > 0 ? 1 : -1);
- z += (z4 > 0 ? 1 : -1);
- /* See if the box is uniformly one side or another. */
- if(z != 4 && z != -4) {
- /* It isn't. Now find the corner of the box with the largest z value that
- has already been plotted, and use the same style used for that node. */
- k = -1000;
- x = -32768;
- if (nodes[j].z > x && nodes[j].style_used !=-1000) {
- k = nodes[j].style_used;
- x = nodes[j].z;
- };
- if (nodes[j+1].z > x && nodes[j+1].style_used !=-1000) {
- k = nodes[j+1].style_used;
- x = nodes[j+1].z;
- };
- if (nodes[j+row_offset+1].z > x && nodes[j+row_offset+1].style_used !=-1000) {
- k = nodes[j+row_offset+1].style_used;
- x = nodes[j+row_offset+1].z;
- };
- if (nodes[j+row_offset].z > x && nodes[j+row_offset].style_used !=-1000) {
- k = nodes[j+row_offset].style_used;
- x = nodes[j+row_offset].z;
- };
- if( k != -1000){
- z = 0; /* To defeat the logic to come. */
- current_style = k;
- (*t->linetype)(current_style);
- };
- };
- /* If k == -1000 then no corner found. I guess it does not matter. */
- };
- if(z > 0 && current_style != plot_info[nplot].above_color) {
- current_style = plot_info[nplot].above_color;
- (*t->linetype)(current_style);
- };
- if(z < 0 && current_style != plot_info[nplot].below_color) {
- current_style = plot_info[nplot].below_color;
- (*t->linetype)(current_style);
- };
- xmin_hl = (sizeof(xleft) == 4 ? 0x7fffffff : 0x7fff );
- xmax_hl = 0;
- clip_move(nodes[j].x,nodes[j].y);
- clip_vector(nodes[j+1].x,nodes[j+1].y);
- clip_vector(nodes[j+row_offset+1].x,nodes[j+row_offset+1].y);
- clip_vector(nodes[j+row_offset].x,nodes[j+row_offset].y);
- clip_vector(nodes[j].x,nodes[j].y);
- nodes[j].style_used = current_style;
- nodes[j+1].style_used = current_style;
- nodes[j+row_offset+1].style_used = current_style;
- nodes[j+row_offset].style_used = current_style;
- MAYBE_LINEPOINT(j);
- MAYBE_LINEPOINT(j+1);
- MAYBE_LINEPOINT(j+row_offset+1);
- MAYBE_LINEPOINT(j+row_offset);
- if( xmin_hl < 0 || xmax_hl > XREDUCE(xright)-XREDUCE(xleft))
- graph_error("Logic error #3 in hidden line");
- /* now mark the area as being filled in the bitmap. These coordinates
- have already been reduced. */
- if (xmin_hl < xmax_hl)
- for(j=xmin_hl;j<=xmax_hl;j++) {
- if (ymin_hl[j] == 0x7fff)
- graph_error("Logic error #2 in hidden line");
- if(pnt[j] == 0) {
- pnt[j] = (unsigned short int *) alloc((unsigned long)y_malloc,"hidden");
- memset(pnt[j], 0, y_malloc);
- };
- if(ymin_hl[j] < 0 || ymax_hl[j] > YREDUCE(ytop)-YREDUCE(ybot))
- graph_error("Logic error #1 in hidden line");
- /* this shift is wordsize dependent */
- indx1 = ymin_hl[j] >> 4;
- indx2 = ymax_hl[j] >> 4;
- mask1 = 0xffff << (ymin_hl[j] & 0xf);
- mask2 = 0xffff >> (0xf-(ymax_hl[j] & 0xf));
- cpnt = pnt[j] + indx1;
- if(indx1 == indx2){
- *cpnt |= (mask1 & mask2);
- } else {
- *cpnt++ |= mask1;
- k = indx1+1;
- while (k != indx2) {
- *cpnt++ = 0xffffU;
- k++;
- };
- *cpnt |= mask2;
- };
- ymin_hl[j]=0x7fff;
- ymax_hl[j]=0;
- };
- };
- };
- free(nodes);
- free(boxlist);
- free(plot_info);
- }
-
-